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Abstract. An isogeny between elliptic curves is an algebraic morphism 
which is a group homomorphism. Many applications in cryptography re- 
quire evaluating large degree isogenies between elliptic curves efficiently. 
For ordinary curves of the same endomorphism ring, the previous best 
known algorithm has a worst case running time which is exponential 
in the length of the input. In this paper we show this problem can be 
solved in subexponential time under reasonable heuristics. Our approach 
is based on factoring the ideal corresponding to the kernel of the isogeny, 
modulo principal ideals, into a product of smaller prime ideals for which 
the isogenies can be computed directly. Combined with previous work of 
Bostan et al., our algorithm yields equations for large degree isogenies 
in quasi-optimal time given only the starting curve and the kernel. 



1 Introduction 

A well known theorem of Tate [29] states that two elliptic curves defined over 
the same finite field ¥ q are isogenous (i.e. admit an isogeny between them) if and 
only if they have the same number of points over ¥ q . Using fast point counting 
algorithms such as Schoof's algorithm and others [9, 25], it is very easy to check 
whether this condition holds, and thus whether or not the curves are isogenous. 
However, constructing the actual isogeny itself is believed to be a hard problem 
due to the nonconstructive nature of Tate's theorem. Indeed, given an ordinary 
curve E/¥ q and an ideal of norm n in the endomorphism ring, the fastest previ- 
ously known algorithm for constructing the unique (up to isomorphism) isogeny 
having this ideal as kernel has a running time of 0(n 3+£ ), except in a certain 
very small number of special cases [4,16,17]. In this paper, we present a new 
probabilistic algorithm for evaluating such isogenies, which in the vast majority 
of cases runs (heuristically) in subexponential time. Specifically, we show that 
for ordinary curves, one can evaluate isogenies of degree n between curves of 
nearly equal endomorphism ring over ¥ q in time less than L q (^, -^) log(n), pro- 
vided n has no large prime divisors in common with the endomorphism ring 
discriminant. Although this running time is not polynomial in the input length, 
our algorithm is still much faster than the (exponential) previous best known 
algorithm, and in practice allows for the evaluation of isogenies of cryptographi- 
cally sized degrees, some examples of which we present here. We emphasize that, 
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in contrast with the previous results of Broker et al. [4], our algorithm is not 
limited to special curves such as pairing friendly curves with small discriminant. 

If an explicit equation for the isogeny as a rational function is desired, our 
approach in combination with the algorithm of Bostan et al. [3] can produce the 
equation in time 0(n 1+e ) given E and an ideal of norm n, which is quasi-optimal 
in the sense that (up to log factors) it is equal to the size of the output. To 
our knowledge, this method is the only known algorithm for computing rational 
function expressions of large degree isogenies in quasi-optimal time in the general 
case, given only the starting curve and the kernel. 

Apart from playing a central role in the implementation of the point count- 
ing algorithms mentioned above, isogenies have been used in cryptography to 
transfer the discrete logarithm problem from one elliptic curve to another [9, 
16,17,20,23,30]. In many of these applications, our algorithm cannot be used 
directly, since in cryptography one is usually given two isogenous curves, rather 
than one curve together with the isogeny degree. However, earlier results [16, 
17,20] have shown that the problem of computing isogenies between a given 
pair of curves can be reduced to the problem of computing isogenies of prime 
degree starting from a given curve. It is therefore likely that the previous best 
isogeny construction algorithms in the cryptographic setting can be improved or 
extended in light of the work that we present here. 

2 Background 

Let E and E' be elliptic curves defined over a finite field ¥ q of characteristic p. An 
isogeny <p ■ E — > E' defined over ¥ q is a non-constant rational map defined over 
¥ q which is also a group homomorphism from E(¥ q ) to E'(¥ q ). This definition 
differs slightly from the standard definition in that it excludes constant maps [27, 
§111.4]. The degree of an isogeny is its degree as a rational map, and an isogeny 
of degree t is called an f-isogeny. Every isogeny of degree greater than 1 can be 
factored into a composition of isogenies of prime degree defined over ¥ q [11]. 

For any elliptic curve E: y 2 + a\xy + a^y = x 3 + a^x 2 + a^x + a@ defined 
over ¥ q , the Frobenius endomorphism is the isogeny ir q : E — »■ £ of degree q 
given by the equation n q (x 7 y) — {x q 1 y q ). The characteristic polynomial of ir q is 
X 2 -tX + q where t = q + l- #E(¥ q ) is the trace of E. 

An endomorphism of E is an isogeny E E defined over the algebraic 
closure ¥ q of ¥ q . The set of endomorphisms of E together with the zero map 
forms a ring under the operations of pointwisc addition and composition; this 
ring is called the endomorphism ring of E and denoted End(-E). The ring End(_E) 
is isomorphic either to an order in a quaternion algebra or to an order in an 
imaginary quadratic field [27, V.3.1]; in the first case we say E is supersingular 
and in the second case we say E is ordinary. 

Two elliptic curves E and E' defined over ¥ q are said to be isogenous over ¥ q if 
there exists an isogeny <fr ■ E — > E' defined over ¥ q . A theorem of Tate states that 
two curves E and E' are isogenous over F g if and only if #E(¥ q ) = #E'(W q ) [29, 
§3]. Since every isogeny has a dual isogeny [27, III. 6.1], the property of being 
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isogenous over ¥ q is an equivalence relation on the finite set of F^-isomorphism 
classes of elliptic curves defined over ¥ q . Moreover, isomorphisms between elliptic 
curves can be classified completely and computed efficiently in all cases [16]. 
Accordingly, we define an isogeny class to be an equivalence class of elliptic 
curves, taken up to F g -isomorphism, under this equivalence relation. 

Curves in the same isogeny class are either all supersingular or all ordinary. 
The vast majority of curves are ordinary, and indeed the number of isomorphism 
classes of supersingular curves is finite for each characteristic. Also, ordinary 
curves form the majority of the curves of interest in applications such as cryp- 
tography. Hence, we assume for the remainder of this paper that we are in the 
ordinary case. 

Let K denote the imaginary quadratic field containing End(£?) , with maximal 
order Ok- For any order O C Ok, the conductor of O is defined to be the 
integer [O k ■ O}. The field K is called the CM field of E. We write c E for 
the conductor of End(E) and for the conductor of T,[ir q ]. It follows from [12, 
§7] that End(-E) = Z + c e Ok and A = c 2 E A K , where A (respectively, A K ) is 
the discriminant of the imaginary quadratic order End(E) (respectively, Ok)- 
Furthermore, the characteristic polynomial has discriminant A v = t 2 — Aq = 
disc(Z[7r 9 ]) = clA K , with c^^c E - [End(£) : Z[tt ? ]]. 

Following [14] and [16], we say that an isogeny </>: E — > E' of prime degree 
£ defined over F 9 is "down" if [End(£) : End(S')] = £, "up" if [End(-E') : 
End(_E)] = I, and "horizontal" if End(_E) = End(E). Two curves in an isogeny 
class are said to "have the same level" if their endomorphism rings are equal. 
Within each isogeny class, the property of having the same level is an equivalence 
relation. A horizontal isogeny always goes between two curves of the same level; 
likewise, an up isogeny enlarges the endomorphism ring and a down isogeny 
reduces it. Since there are fewer elliptic curves at higher levels than at lower 
levels, the collection of elliptic curves in an isogeny class visually resembles a 
"pyramid" or a "volcano" [14], with up isogenies ascending the structure and 
down isogenies descending. If we restrict to the graph of ^-isogenies for a single 
£, then in general the i?-isogeny graph is disconnected, having one ^-volcano 
for each intermediate order Z[n q ] C O C Ok such that O is maximal at I 
(meaning £ j [Ok ■ O]). The "top level" of the class consists of curves E with 
End(_E) = Ok, and the "bottom level" consists of curves with End(E) = Z[n q ]. 

We say that £ is an Elkies prime [2, p. 119] if £ \ ce and {—) 7^ — 1> or 
cquivalently if and only if E admits a horizontal isogeny of degree £. The number 
of ^-isogenics of each type can easily be determined explicitly [14,16,21]. In 
particular, for all but the finitely many primes £ dividing [Ok ■ Z[7r g ]], we have 
that every rational ^-isogeny admitted by E is horizontal. 

3 The Broker-Charles-Lauter algorithm 

Our algorithm is an extension of the algorithm developed by Broker, Charles, 
and Lauter [4] to evaluate large degree isogenies over ordinary elliptic curves with 
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endomorphism rings of small class number, such as pairing- friendly curves [15]. 
In this section we provide a summary of their results. 

The following notation corresponds to that of [4]. Let E/¥ q be an ordi- 
nary elliptic curve with endomorphism ring End(£?) isomorphic to an imaginary 
quadratic order Oa of discriminant A < 0. Identify End(-E) with Oa via the 
unique isomorphism t such that i* (x)lo = xuj for all invariant differentials to and 
all x <G Oa- Then every horizontal separable isogeny on E of prime degree £ 
corresponds (up to isomorphism) to a unique prime ideal £ C O a of norm £ for 
some Elkies prime £. We denote the kernel of this isogeny by E[£\. Any two dis- 
tinct isomorphic horizontal isogenies (i.e., pairs of isogenies where one is equal to 
the composition of the other with an isomorphism) induce different maps on the 
space of differentials of E, and a separable isogeny is uniquely determined by the 
combination of its kernel and the induced map on the space of differentials. A 
normalized isogeny is an isogeny <p: E — > E' for which <ff(u>E') — where loe 
denotes the invariant differential of E. Algorithm 1 (identical to Algorithm 4.1 
in [4]) evaluates, up to automorphisms of E, the unique normalized horizontal 
isogeny of degree £ corresponding to a given kernel ideal £ C Oa- 

The following theorem, taken verbatim from [4] , shows that the running time 
of Algorithm 1 is polynomial in the quantities \og(£), \og(q), n, and \A\. 

Theorem 3.1. Let E/¥ q be an ordinary elliptic curve with Frobenius n q , given 
by a Weierstrass equation, and let P e E(F q ™) be a point on E. Let A = 
disc(End(_E)) be given. Assume that [End(_E) : Z[7r g ]] and #E(¥ q n) are coprime, 
and let £ = (£, c + dir q ) be an End(£ l ) -ideal of prime norm £ ^ char(F g ) not 
dividing the index [End(.E') : Z[V g ]]. Algorithm 1 computes the unique elliptic 
curve E' such that there exists a normalized isogeny (f>: E — > E' with kernel 
E[£\. Furthermore, it computes the x-coordinate of <p(P) if End(_E) does not 
equal Z[i] or 1\C,z] and the square, respectively cube, of the x-coordinate of <f>(P) 
otherwise. The running time of the algorithm is polynomial in log(^), \og(q), n 
and \A\. 

4 A subexponential algorithm for evaluating horizontal 
isogenies 

As was shown in Sections 2 and 3, any horizontal isogeny can be expressed as a 
composition of prime degree isogenies, one for each prime factor of the kernel, 
and any prime degree isogeny is a composition of a normalized isogeny and 
an isomorphism. Therefore, to evaluate a horizontal isogeny given its kernel, it 
suffices to treat the case of horizontal normalized prime degree isogenies. 

Our objective is to evaluate the unique horizontal normalized isogeny on a 
given elliptic curve E/¥ q whose kernel ideal is given as £ = (£, c+dn q ), at a given 
point P e E(¥ q n), where £ is an Elkies prime. As in [4], we must also impose the 
additional restriction that £ \ [End(_E) : Z[7r q ]]; for Elkies primes, an equivalent 
restriction is that I \ [Ok ■ %[xq]], but we retain the original formulation for 
consistency with [4] . 
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Algorithm 1 The Broker-Charles-Lauter algorithm 

Input: A discriminant A, an elliptic curve E/¥ q with End(i5) = Oa and a point P £ 

E{¥ q n) such that [End(_E) : Z[iv q ]] and #E(W q ™) are coprime, and an End(_E)-ideal 

£ = (£, c+dir q ) of prime norm £ 7^ char(F q ) not dividing the index [End(_E) : Z[7r 9 ]]. 
Output: The unique elliptic curve E' admitting a normalized isogeny cf>: E — > E' with 

kernel E[£], and the ^-coordinate of 4>(P) for A —3,-4 and the square (resp. 

cube) of the ^-coordinate otherwise. 
1: Compute the direct sum decomposition Pic(Ozi) = &>{[Ii\) of Pic(Oa) into cyclic 

groups generated by the degree 1 prime ideals h of smallest norm that are coprime 

to the product p ■ #E(¥ q ^) • [End(_E) : Z[tt,]]. 
2: Using brute force , find ei, e2, . . . , such that [£] = [II 1 ] ■ [I^ 2 ] • • ■ [Il k ]- 
3: Find a (using Cornacchia's algorithm) and express £ = I^ 1 ■ J| 2 • • • ■ (a). 
4: Compute a sequence of isogenies (0i, . . . , cf> s ) such that the composition <p c : E — > 

E c has kernel E[I^ ■ J 2 e2 • • • I e k k ] using the method of [4, § 3]. 
5: Evaluate cj> c (P) e E c (W q n). 

6: Write a = (u + vir q ) / (zm) . Compute the isomorphism 77: E c ^ E' with rf{oj E ') = 

(u/ zm)ujE a - Compute Q = ri(<p c (P)). 
7: Compute (zm) -1 mod #E(F q n), and compute R = ((zm) _1 (it + viv q ))(Q). 
8: Put r = x(Ry° Ar/2 and return (E',r). 



In practice, one is typically given £ instead of £, but since it is easy to cal- 
culate the list of (at most two) possible primes £ lying over £ (cf. [6]), these two 
interpretations are for all practical purposes equivalent, and we switch freely 
between them when convenient. When £ is small, one can use modular polyno- 
mial based techniques [4, §3.1], which have running time 0(£ 3 \og(£ ) 4+£ ) [13]. 
However, for isogeny degrees of cryptographic size (e.g. 2 160 ), this approach 
is impractical. The Broker-Charles-Lauter algorithm sidesteps this problem, by 
using an alternative factorization of £. However, the running time of Broker- 
Charles-Lauter is polynomial in |Z\|, and therefore even this method only works 
for small values of \A\. In this section we present a modified version of the 
Broker-Charles-Lauter algorithm which is suitable for large values of \A\. 

We begin by giving an overview of our approach. In order to handle large 
values of there are two main problems to overcome. One problem is that we 
need a fast way to produce a factorization 

Z = IZ 1 I?--I?-(a) (1) 

as in lines 2 and 3 of Algorithm 1. The other problem is that the exponents 
in Equation (1) need to be kept small, since the running times of lines 3 and 4 
of Algorithm 1 arc proportional to J^i \ e i\ Norm(/i) 2 . The first problem, that 
of finding a factorization of £, can be solved in subexponential time using the 
index calculus algorithm of Hafner and McCurley [18] (see also [6, Chap. 11]). 

1 Broker, Charles, and Lauter mention that this computation can be done in "various 
ways" [4, p. 107], but the only explicit method given in [4] is brute force. The use 
of brute force limits the algorithm to elliptic curves for which \A\ is small, such as 
pairing-friendly curves. 
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Algorithm 2 Computing a factor base 
Input: A discriminant A, a bound N. 

Output: The set I consisting of split prime ideals of norm less than N, together with 
the corresponding set 7 of quadratic forms. 
1: Set 7i- 0. 
2: SetX^0. 

3: Find all primes p < N such that ( j) = 1. Call this set P. Let k = \P\. 

4: For each prime Pi £ P, find an ideal pi of norm pi (using Cornacchia's algorithm). 

5: For each i, find a quadratic form fi — [(pi,bi,Ci)] corresponding to pi in C\(Oa), 

using the technique of [26, §3]. 
6: Output 1 = {pi,p 2 , • • • ,pk} and 7 = {/i,/ 2 , ■■■,fk}- 



To resolve the second problem, we turn to an idea which was first introduced by 
Galbraith et. al [17], and recently further refined by Bisson and Sutherland [1]. 
The idea is that, in the process of sieving for smooth norms, one can arbitrar- 
ily restrict the input exponent vectors to sparse vectors (ei, e<i, e^) such that 
J2i \ e i\N{Ii) 2 is kept small. This restriction is implemented in line 6 of Algo- 
rithm 3. As in [1], one then assumes heuristically that the imposition of this 
restriction does not affect the eventual probability of obtaining a smooth norm 
in the Hafner and McCurley algorithm. Note that, unlike the input exponents, 
the exponents appearing in the factorizations of the ensuing smooth norms (that 
is, the values of yi in Algorithm 3) are always small, since the norm in question 
is derived from a reduced quadratic form. 

We now describe the individual components of our algorithm in detail. 

4.1 Finding a factor base 

Let C\(Oa) denote the ideal class group of O a- Algorithm 2 produces a factor 
base consisting of split primes in O a of norm less than some bound N. The 
optimal value of N will be determined in Section 4.4. 

4.2 "Factoring" large prime degree ideals 

Algorithm 3, based on the algorithm of Hafner and McCurley, takes as input a 
discriminant A, a curve E, a prime ideal £ of prime norm I in Oa, a smoothness 
bound N, and an extension degree n. It outputs a factorization 

£ = /- If ■■■II-- (a) 

as in Equation 1, where the Ij's are as in Algorithm 1, the exponents ei are 
positive, sparse, and small (i.e., polynomial in N), and the ideal (a) is a principal 
fractional ideal generated by a. 

4.3 Algorithm for evaluating prime degree isogenies 

The overall algorithm for evaluating prime degree isogenies is given in Algo- 
rithm 4. This algorithm is identical to Algorithm 1, except that the factoriza- 
tion of £ is performed using Algorithm 3. To maintain consistency with [4], we 
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Algorithm 3 "Factoring" a prime ideal 

Input: A discriminant A, an elliptic curve E/¥ q with End(iS) = Oa, a smoothness 

bound N, a prime ideal £ of norm £ in Oa, an extension degree n. 
Output: Relation of the form £ = (a) • n i=1 where (a) is a fractional ideal, Ii 
are as in Algorithm 1, and e, > are small and sparse. 
1: Run Algorithm 2 on input A and JV to obtain X = {pi,p2, ■ ■ ■ ,pk} and J 7 = 

{/i, /2, • • • , /fe}- Discard any primes dividing p ■ #£(F 9 n) • [End(_B) : Z[n q ]]. 
2: Set pi <— Norm(p»). (These values are also calculated in Algorithm 2.) 
3: Obtain the reduced quadratic form [£] corresponding to the ideal class of £. 
4: repeat 

5: for i = 1, . . . , k do 

6: Pick exponents Xi in the range [0, (N/pi) 2 ] such that at most ko are nonzero, 

where ko is a global absolute constant (in practice, fco = 3 suffices) . 
7: end for 

8: Compute the reduced quadratic form o = (a, b, c) for which the ideal class [a] is 

equivalent to [£] • n»=i fi* ■ 
9: until The integer a factors completely into the primes pi, and the relation derived 

from [a] = [£] ■ fTjli fi 1 contains fewer than y/log(\A\/3)/z nonzero exponents. 
10: Write a = U^iPT- 
11: for i=l, . . . , k do 

12: Using the technique of Seysen ([26, Theorem 3.1]), determine the signs of the 

exponents j/; = ±m for which o = n»=i fP ■ 
13: Let e, = j/» — a;<. (These exponents satisfy [£] = n»=i fi' ■) 
14: if a > then 
15: Set Ji <- pi 

16: else 
17: Set Ii «- pi 

18: end if 
19: end for 

20: Compute the principal ideal I = £ • F]Li ^i^'- 

21: Using Cornacchia's algorithm, find a generator j3 G of I. 

22: Sctm^ntiPi ed and q 4 ^ . 

23: Output £=(«)• /j ei1 ■ 4 e2 ' • • • I l k ekl . 



have included the quantities A and End(E') as part of the input to the algo- 
rithm. However, we remark that these quantities can be computed from E/¥ q 
in L 9 (|, ^) operations using the algorithm of Bisson and Sutherland [1], even 
if they are not provided as input. 

4.4 Running time analysis 

In this section, we determine the theoretical running time of Algorithm 4, as well 
as the optimal value of the smoothness bound N to use in line 1 of the algorithm. 
As is typical for subexponential time factorization algorithms involving a factor 
base, these two quantities depend on each other, and hence both are calculated 
simultaneously. 
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Algorithm 4 Evaluating prime degree isogenies 

Input: A discriminant A, an elliptic curve E/¥ q with End(i5) = Oa and a point P £ 
E{¥ q n) such that [End(_E) : Z[-ir q ]] and #_E(F 9 n) are coprime, and an End(_E)-ideal 
£ = (I, c+d-Kq) of prime norm £ 7^ char(F q ) not dividing the index [End(_E) : Z[7r 9 ]]. 

Output: The unique elliptic curve E' admitting a normalized isogeny cf>: E — > E' with 
kernel E[£], and the ^-coordinate of 4>(P) for A —3,-4 and the square (resp. 
cube) of the ^-coordinate otherwise. 
1: Choose a smoothness bound N (see Section 4.4). 

2: Using Algorithm 3 on input (A, E, N, £, n), obtain a factorization of the form 

£ = /r i -i5 2 "-c-(«)- 

3: Compute a sequence of isogenies (</>i, . . . , S ) such that the composition <^> c : E — > 

£ c has kernel S^ 1 • J 2 62 • • • I e k k ] using the method of [4, § 3]. 
4: Evaluate cj> c (P) G £ c (F 9 n). 

5: Write a = (u + vn q )/(zm). Compute the isomorphism q: E c — > E' with t]*(u}e') = 

(u/ zm)ujE a - Compute Q = r)(<f) c (P)). 
6: Compute (zm)^ 1 mod#_E(F q n), and compute i? = ((m) _1 (ii + «7r 9 ))(Q). 
7: Put r = ^(i?) 10 ^ 1 */ 2 and return (E',r). 



As in [9], we define 2 L n (a,c) by 

i„( a , C )=0(exp((c + o(l))(log(n)) Q (log(log(n))) 1 - Q )). 

The quantity L n (a,c) interpolates between polynomial and exponential size as 
a ranges from to f . We set N = (|, z) for an unspecified value of z, and in 
the following paragraphs we determine the optimal value of z which minimizes 
the running time of Algorithm 4. (The fact that a — | is optimal is clear from 
the below analysis, as well as from prior experience with integer factorization 
algorithms.) For convenience, we will abbreviate (a, c) to L(a, c) throughout. 

Line 2 of Algorithm 4 involves running Algorithm 3, which in turn calls 
Algorithm 2. As it turns out, Algorithm 2 is almost the same as Algorithm 11.1 
from [6], which requires L(^, z) time, as shown in [6]. The only difference is that 
we add an additional step where we obtain the quadratic form corresponding to 
each prime ideal in the factor base. This extra step requires 0(log(Norm(/)) 1+e ) 
time for a prime ideal /, using Cornacchia's Algorithm [19]. Thus, the overall 
running time for Algorithm 2 is bounded above by 

L{\,z) ■ log(L(i, z)) 1+e = i(i, z). 

Line 2 of Algorithm 3 takes log(£) time using standard algorithms [12]. The 
loop in lines 4-9 of Algorithm 3 is very similar to the FindRelation algorithm 
in [1], except that we only use one discriminant, and we omit the requirement 
that #R/Di > #R/D 2 (which in any case is meaningless when there is only 
one discriminant). Needless to say, this change can only speed up the algorithm. 
Taking fi = V2 z in [1, Prop. 6], we find that the (heuristic) expected running 
time of the loop in lines 4-9 of Algorithm 3 is L(^,j-). 

2 The definition of L n (a, c) in [6] differs from that of [9] in the o(l) term. We account 
for this discrepancy in our text. 
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The next step in Algorithm 3 having nontrivial running time is the computa- 
tion of the ideal product in line 20. To exponentiate an element of an arbitrary 
semigroup to a power e requires 0(log e) semigroup multiplication operations [10, 
§1.2]. To multiply two ideals I and J in an imaginary quadratic order (via com- 
position of quadratic forms) requires 0(max(log(Norm(/)), log(Norm( J))) 1+e ) 
bit operations using fast multiplication [24, §6]. Each of the expressions |/i|' ei ' 
therefore requires 0(log|ej|) ideal multiplication operations to compute, with 
each individual multiplication requiring 



0((| ei | log(Norm(/ l ))) 1+e ) = O - logfa) = 0(7V 2+£ ) 




bit operations, for a total running time of (log ei)0(N 2+e ) — L(^,2z) for each i. 
This calculation must be performed once for each nonzero exponent e,. By 
line 9, the number of nonzero exponents appearing in the relation is at most 
-\/log(|Z\|/3)/z, so the amount of time required to compute all of the |/»|' e< ' for 
all i is (y/log(\A\/3)/z)L(±,2z) = L(\,2z). Afterward, the values \I,\^ must 
all be multiplied together, a calculation which entails at most ^/log(|Z\|/3)/z 
ideal multiplications where the log- norms of the input multiplicands are bounded 
above by 

logNorm(/] ei1 ) = |e»| logNorm(/ 4 ) < (^j log Pi < N 2 = L(^,2z), 

and thus each of the (at most) y/log(\A\/3)/z multiplications in the ensuing 
product can be completed in time at most {^\og{\A\/3) / z)L{^, 2z) — L(^,2z). 
Finally, we must multiply this end result by £, an operation which requires 
0(max(log^, L(|, 2z)) 1+£ ) time. All together, the running time of step 20 is 
L(i,2z) + 0(max(log£,L(i,2z)) 1+£ ) = max((log^) 1+e , L{\, 2z)), and the norm 
of the resulting ideal I is bounded above by I ■ exp(L(|, 2z)). 

Obtaining the generator /? of J in line 21 of Algorithm 3 using Cornacchia's 
algorithm requires 



0(log(Norm(/)) 1+£ ) = {\ogl + L{\, 2z)) 



l+e 



time. We remark that finding (3 given I is substantially easier than the usual Cor- 
nacchia's algorithm, which entails finding (i given only Norm(J). The usual algo- 
rithm requires finding all the square roots of A modulo Norm(J), which is very 
slow when Norm(I) has a large number of prime divisors. This time-consuming 
step is unnecessary when the ideal I itself is given, since the embedding of the 
ideal I in End(-E') already provides (up to sign) the correct square root of A 
mod /. A detailed description of this portion of Cornacchia's algorithm in the 
context of the full algorithm, together with running time figures specific to each 
sub-step, is given by Hardy et al. [19]; for our purposes, the running time of a 
single iteration of Step 6 in [19, §4] is the relevant figure. This concludes our 
analysis of Algorithm 3. 
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Returning to Algorithm 4, we find that (as in [4]) the computation of the indi- 
vidual isogenics (j>i in line 3 of Algorithm 4 is limited by the time required to com- 
pute the modular polynomials 4> n (x,y). Using the Chinese remainder theorem- 
based method of Broker et al. [5] , these polynomials can be computed mod q in 
time 0(n 3 log 3+e (n)), and the resulting polynomials require 0(n 2 (log 2 ro + logq)) 
space. For each ideal ij, the corresponding modular polynomial of level pi only 
needs to be computed once, but the polynomial once computed must be eval- 
uated, differentiated, and otherwise manipulated a times, at a cost of 0(p 2+£ ) 
field operations in ¥ q per manipulation, or 0(p 2+e )(\og q) 1+£ bit operations using 
fast multiplication. The total running time of line 3 is therefore 



< 0(N 3+£ ) + vMgEg) jv2+e (logg) i+ £ = L( i ; 3z) + L{ i ; 2z)(\o gq y+t. 
Similarly, the evaluation of <p c in line 4 requires 



^\e t \p^ = L(l2z) 



field operations in F q n, which corresponds to L(^,2z)(\ogq n ) 1+£ bit operations 
using fast multiplication. 

Combining all the above quantities, we obtain a total running time of 



(algorithm 2) 

+ L{\,±) (lines 4-9, algorithm 3) 

+ max((log£) 1 + £ ,L(±,2;j)) (line 20, algorithm 3) 

+ (log^ + L(±,2z)) 1+£ (line 21, algorithm 3) 

+ L(i,3z) + L(i,2z)(logq) 1+e (line 3, algorithm 4) 

+ L{\, 2z){\ogq n ) 1+£ (line 4, algorithm 4) 

= L(|, + (log^ + L{\, 2z)f+ £ + L&3z) + L{\, 2z){\ogq-f+ £ . 

When \A\ is large, we may impose the reasonable assumption that log(^) <C 
L(i, z) and log(g") <C L(^, z). In this case, the running time of Algorithm 4 is 
dominated by the expression L(^, j^) + L(^,3z), which attains a minimum at 
z = Taking this value of z, we find that the running time of Algorithm 4 

is equal to ^). Since the maximum value of \A\ < \A n \ = Aq — t 2 is Aq, 

we can alternatively express this running time as simply L q (^, -^). 

In the general case, log(^) and log(g n ) might be non-negligible compared to 
L{\,z). This can happen in one of two ways: either \A\ is small, or (less likely) 
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£ is very large and/or n is large. When this happens, we can still bound the 
running time of Algorithm 4 by taking z = 7^ in the foregoing calculation, 
although such a choice may fail to be optimal. We then find that the running 
time of Algorithm 4 is bounded above by 

(logW + L{\, + L(§, % + L(±, ^)(log g«)^. 

We summarize our results in the following theorem. 

Theorem 4.1. Let E/¥ q be an ordinary elliptic curve with Frobenius ir q , given 
by a Weierstrass equation, and let P G E(¥ q n) be a point on E. Let A = 
disc(End(_E)) be given. Assume that [End(_E) : Z[7r g ]] and #E(F qn ) are coprime, 
and let £ = (£,c+ dn q ) be an F,nd(E)-ideal of prime norm £ ^ char(F g ) not 
dividing the index [End(.E') : Z[iT q ]]. Under the heuristics of [1, $4-]i Algorithm 4 
computes the unique elliptic curve E' such that there exists a normalized isogeny 
<p: E — > E' with kernel E[£\. Furthermore, it computes the x-coordinate of <j>(P) 
if End(_E) does not equal Z[i] or ^[£3] and the square, respectively cube, of the 
x-coordinate of <f>(P) otherwise. The running time of the algorithm is bounded 
above by 

(\og(£) + L{\, ^))^ + L(|, % + L(|, ^)(log q-)^. 

The running time of the algorithm is subexponential in log \A\, and polynomial 
in log(£), log(g), and n. 



5 Examples 

5.1 Small example 

Let p = 10 10 + 19 and let E/¥ p be the curve y 2 = a; 3 + 15.T + 129. Then E(¥ p ) has 
cardinality 10000036491 = 3 • 3333345497 and trace t = -36471. To avoid any 
bias in the selection of the prime I , we set £ to be the smallest Elkies prime of 
E larger than p/2, namely £ — 5000000029. We will evaluate the x-coordinate of 
4>{P), where (f) is an isogeny of degree £, and P is chosen arbitrarily to be the point 
(5940782169,2162385016) e E(¥ p ). We remark that, although this example is 
designed to be artificially small for illustration purposes, the evaluation of this 
isogeny would already be infeasible if we were using prior techniques based on 
modular functions of level £. 

The discriminant A of E is A = t 2 - Ap = -38669866235. Set w = ±±y^2 
and O = Oa- The quadratic form (5000000029,-2326859861,270713841) rep- 
resents a prime ideal £ of norm £, and we show how to calculate the isogeny <j> 
having kernel corresponding to E[£\. Using an implementation of Algorithm 3 
in MAGMA [22], we find immediately the relation £ = (£) • p ig ■ p§f where 
P = 588048307603210005w - 235788727470005542279904, m = 19 • 31 24 , p i9 = 
(19, 2w + 7), and p 3i = (31, 2w + 5). Using this factorization, we can then eval- 
uate (f>: E — > E' using the latter portion of Algorithm 4. We find that E' is 
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the curve with Weierstrass equation y 2 = x 3 + 3565469415.x + 7170659769, and 
4>{P) = (7889337683, ±3662693258). We omit the details of these steps, since 
this portion of the algorithm is identical to the algorithm of Broker, Charles 
and Lauter, and the necessary steps are already extensively detailed in their 
article [4]. 

We can check our computations for consistency by performing a second com- 
putation, starting from the curve E' : y 2 = x 3 + 3565469415x+ 7170659769, the 
point P' = (7889337683, 3662693258) € E'(¥ p ), and the conjugate ideal £, which 
is represented by the quadratic form (5000000029, 2326859861, 270713_841). Let 
<p: E' — > E" denote the unique normalized isogeny with kernel £"[£]. Up to 
a normalization isomorphism t: E — > E", the isogeny <f> should equal the dual 
isogeny (f> of 0, and the composition </>(0(P)) should yield l(£P). Indeed, upon 
performing the computation, we find that E" has equation 



which is isomorphic to E via the isomorphism i: E — > E" defined by t(x,y) = 
(x/£ 2 ,y/£ 3 ), and 

4>{4>{P)) = (3163843645,8210361642) = (5551543736/^ 2 , 6305164567/f 5 ), 
in agreement with the value oilP, which is (5551543736,6305164567). 

5.2 Medium example 

Let E be the ECCp-109 curve [8] from the Certicom ECC Challenge [7], with 
equation y 2 = x 3 + ax + b over ¥ p where 



As before, to avoid any bias in the choice of £, we set I to be the least Elkies prime 
greater than p/2, and we define w = where A = disc(End(S)). Let £ be 

the prime ideal of norm £ in End(-E) corresponding to the reduced quadratic form 
{£,b,c) of discriminant A, where b = -105137660734123120905310489472471. 
For each Elkies prime p, let p p denote the unique prime ideal corresponding to 
the reduced quadratic form (p, b, c) where b > 0. Our smoothness bound in this 
case is N — L(i, ^^) w 200. Using Sutherland's smoothrelation package [28], 
which implements the FindRelation algorithm of [1], one finds in a few seconds 
(using an initial seed of 0) the relation £ 3, where 



y 2 = x 3 + (lb/£ 4 )x + (l29/£ 6 ), 



p = 564538252084441556247016902735257 
a = 321094768129147601892514872825668 



b = 430782315140218274262276694323197 




3 - P7 2 P13°P23P47P73P103P179P191 

to = 7 72 13 100 23 14 47 2 73 2 103 1 179 1 191 1 
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and 

/3 = 3383947601020121267815309931891893555677440374614137047492987151\ 
2226041731462264847144426019711849448354422205800884837 
- 1713152334033312180094376774440754045496152167352278262491589014\ 
097167238827239427644476075704890979685 • w 

We find that the codomain E' of the normalized isogeny <p ■ E — > E' of kernel 
E[£] has equation y 2 = x 3 + a'x + b' where 

a' = 84081262962164770032033494307976 
b' = 506928585427238387307510041944828 

and that the base point 

P = (97339010987059066523156133908935, 149670372846169285760682371978898) 

of E given in the Certicom ECC challenge has image 

(450689656718652268803536868496211, ±345608697871189839292674734567941). 

under cf). As with the first example, we checked the computation for consistency 
by using the conjugate ideal. 

5.3 Large example 

Let E be the ECCp-239 curve [8] from the Certicom ECC Challenge [7]. Then 
E has equation y 2 = x 3 + ax + b over F p where 

p= 862591559561497151050143615844796924047865589835498401307522524859467869 
a = 820125117492400602839381236756362453725976037283079104527317913759073622 
6=545482459632327583111433582031095022426858572446976004219654298705912499 

Let £ be the prime ideal whose norm is the least Elkies prime greater than 
p/2 and whose ideal class is represented by the quadratic form (£,b,c) with 
b > 0. We have = L(i, s=y 5000, and one finds in a few hours using 
smoothrelation [28] that £ is equivalent to 

r» —2 2 — 2 — _ _ _ —18 — — — — 6 —5 — 

^ = P7pllpl9p37p7lpl31p211p389p433p467p859p863pl019pll51pl597p2143p2207p3359 

where each ideal p p is represented by the reduced quadratic form (p, b, c) having 
b > (this computation can be reconstructed with [28] using the seed 7). The 
quotient £/3 is generated by /3/ra where m — Norm(J) and ft is 

-923525986803059652225406070265439117913488592374741428959120914067053307\ 
4585317 - 917552768623818156695534742084359293432646189962935478129227909w. 

Given this relation, evaluating isogenies of degree i is a tedious but routine com- 
putation using Elkies- Atkin techniques [4, §3.1]. Although we do not complete 
it here, the computation is well within the reach of present technology; indeed, 
Broker et al. [5] have computed classical modular polynomials mod p of level up 
to 20000, well beyond the largest prime of 3389 appearing in our relation. 



14 



David Jao and Vladimir Soukharev 



6 Related work 

Bisson and Sutherland [1] have developed an algorithm to compute the endo- 
morphism ring of an elliptic curve in subexponential time, using relation-finding 
techniques which largely overlap with ours. Although our main results were ob- 
tained independently, we have incorporated their ideas into our algorithm in 
several places, resulting in a simpler presentation as well as a large speedup 
compared to the original version of our work. 

Given two elliptic curves E and E' over F 9 admitting a normalized isogeny 
<fi: E — > E' of degree £, the equation of as a rational function contains 0(£) 
coefficients. Bostan et al. [3] have published an algorithm which produces this 
equation, given E, E' , and I. Their algorithm has running time 0(£ 1+£ ), which 
is quasi-optimal given the size of the output. Using our algorithm, it is possible 
to compute E' from E and £ in time log(^)Li z \i (_, ^) for large I. Hence the 
combination of the two algorithms can produce the equation of <p within a quasi- 
optimal running time of 0(£ 1+£ ), given only E and £ (or E and £), without the 
need to provide E' in the input. 
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